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ABSTRACT 


Several  location  methods  developed  at  SDAC  are  evaluated  by  compar¬ 
ing  location  errors  of  these  methods  against  errors  computed  with  the 
standard  method.  The  new  methods  include  : 

Using  P  travel  times  in  laterally  heterogeneous  media; 

A  station  travel  time  residual  correlation  matrix  in  the  normal 

equations; 

Location  with  simultaneous  determination  of  Pn,  Pg  and  Lg  veloci¬ 
ties; 

Combination  of  the  correlation  matrix  method  and  the  simultaneous 

determination  method. 

All  the  methods  evaluated  in  this  report  appear  to  give  smaller 
absolute  epicenter  and  depth  errors  than  the  standard  method  of  locating 
events.  Furthermore,  these  methods  do  not  use  station  corrections; 
therefore  the  locations  obtained  are  not  station-  or  region-dependent. 
They  can  be  used  to  locate  seismic  events  in  any  source  region,  with  any 
set  of  stations. 

When  average  location  errors  are  compared  against  errors  with  the 
standard  method,  the  correlation  matrix  method  is  about  3  km  better, 
e.g.,  about  17  kilometers  as  compared  to  20  kilometers;  the  regional 
models  method  is  3  to  6  km  better,  e.g.,  about  5  kilometers  as  compared 
to  9  kilometers,  depending  on  whether  near-regional  stations  are  used. 
The  average  location  errors  from  the  true  epicenters  using  regional  data 
are  approximately  6  km. 

The  addition  of  Lg  arrivals  to  locate  events  did  not  result  in 
better  locations.  Typically,  we  added  6  Lg  arrival  times  to  a  total  of 
15  Pn  and  Pg  arrival  times.  The  lack  of  improvement  may  be  attributed 
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to  a  greater  variance  in  Lg  travel  times. 

When  events  were  located  with  both  teleseismic  P  and  regional 
phases,  the  weighting  factors  for  these  phases  played  an  important  role. 
When  all  phases  were  weighted  equally,  location  errors  were  smaller  than 
when  the  P  phases  were  weighted  more  heavily.  Depth  errors  using  all 
phases  were  similar  to  those  computed  with  teleseismic  P  arrivals  only. 
This  is  to  be  expected  because  the  travel  times  of  regional  phases  are 
affected  weakly,  if  at  all,  by  event  depth. 

The  six  most  thoroughly  analyzed  events  were  FAULTLESS,  RULISON, 
SHOAL,  GASBUGGY,  SALMON,  and  GNOME.  The  median  number  of  P  arrivals  was 
74;  of  Pn.  10;  of  Pg,  1;  and  of  Lg,  4.  The  median  location  error  using 
the  standard  method  was  about  20  kilometers  and  the  median  depth  error 
was  also  about  20  kilometers.  By  adding  the  regional  phase  arrivals  and 
weighting  each  phase  equally,  these  numbers  were  reduced  to  3.6  and  3.7 
kilometers  respectively.  Furthermore,  while  the  true  epicenter  lay  out¬ 
side  the  95Z  confidence  ellipse  5  out  of  6  times  in  the  standard  case, 
it  lay  within  the  ellipse  in  all  cases  in  the  improved  analyses.  These 
results  represent  a  remarkable  improvement  in  absolute  location  accu¬ 
racy.  Most  of  the  improvement  is  due  to  the  use  of  regional  phases. 
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1 .  INTRODUCTION 


For  the  past  several  years,  we  have  been  investigating  various 
methods  for  improving  seismic  location  estimates.  Improvements  can  be 
made  when  the  Earth's  lateral  heterogeneities  are  properly  handled. 
However,  this  is  not  a  simple  task  because  large  fluctuations  in  layer 
velocities  are  observed  over  very  small  distances. 

Seismic  events  are  traditionally  located  at  teleseismic  distances 
with  teleseismic  P  arrivals.  The  use  of  later  phases  is  somewhat  res¬ 
tricted  for  various  reasons.  For  example,  PKP  arrivals  are  often  used 
together  with  P.  but  only  one  branch,  DF,  has  typically  been  used  in 
location  programs.  On  the  other  hand,  the  BC  branch  has  the  greatest 
amplitudes,  and  BC  arrivals  are  often  observed  instead  of  DF  arrivals. 
pP  arrivals  are  very  useful  in  identifying  the  event  depth,  but  observa¬ 
tion  of  pP  is  often  not  easy  because  of  the  presence  of  the  coda.  We 
have  developed  a  program  which  has  two  distinct  options  for  locating 
seismic  events  using  teleseismic  P.  The  first  option  computes  source  to 
station  travel  times  by  using  &  priori  knowledge  of  regional  Earth 
models.  The  second  option  uses  a  correlation  matrix  in  the  normal  equa¬ 
tions  in  which  the  matrix  elements  are  calculated  from  a  simple  function 
of  inter-station  distance. 

Studies  in  recent  years  have  shown  that  the  crust  of  the  Earth  is 
probably  more  heterogeneous  than  we  previously  thought.  Lateral  changes 
in  layer  velocities  can  be  on  the  order  of  0.1  km/sec  instead  of  0.01 
km/ sec.  Under  such  conditions,  it  can  be  very  difficult  to  obtain  accu- 
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rate  location  estimates  using  Pn  and  Pg  arrivals.  Although  some  workers 
have  shown  that  accurate  locations  can  be  obtained  by  using  three- 
dimensional  ray-tracing  programs  to  compute  travel  times,  such  methods 
require  exact  knowledge  of  the  structure  of  the  region.  Thus, these 

methods  are  not  applicable  to  an  area  of  unknown  geologic  structures, 

and  seismic  event  location  in  such  areas  is  often  necessary.  Our  inves¬ 

tigations  in  seismic  location  using  Pn  and  Pg  resulted  in  two  possible 
methods.  Both  begin  the  location  with  a  simple  standard  Earth  model  and 
modify  layer  velocities  as  the  location  proceeds.  The  methods  determine 
layer  velocities  by  either  successive  determinations  (SUCCESS)  or  simul¬ 
taneously  with  seismic  locations  (SIMUL),  and  both  give  much  better 
results  than  the  traditional  location  method. 

These  studies  have  resulted  in  some  useful  options  for  location 
programs.  However,  these  options  had  not  been  collected  in  a  program 
which  could  easily  use  them  all.  Our  objective  was  to  create  a  simpli¬ 
fied  location  program  with  all  the  available  options  combined.  This 

program  was  developed  on  the  VAX  11/780  for  the  Regional  Event  Location 
System  (RELS).  This  report  describes  the  program  developments  in  RELS, 
and  some  results  of  an  evaluation  of  the  various  location  methods 
developed  in  the  previous  years. 
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2.  PROGRAM  DEVELOPMENTS 

Conversion  of  the  location  program  from  the  TS44  system  to  the 
RELS  on  VAX  11/780  seemed  to  be  a  straightforward  task,  but  it  took 
longer  than  expected.  Part  of  the  reason  was  a  lack  of  programming  sup¬ 
port  and  some  delays  in  system  development.  However,  some  significant 
changes  in  the  location  program  were  made  during  this  period.  The  fol¬ 
lowing  highlights  some  of  these  changes. 

2.  ,1_.  Options  in  the  Location  Program 

Four  options  in  the  location  program  were  developed  in  the  past. 
They  are: 

(1)  ACLOC:  Locations  using  models  of  the  laterally  heterogeneous  crust 
and  upper  mantle.  Travel  times  were  not  computed  by  ray  tracing, 
but  by  adding  travel  time  corrections  to  the  Herrin  '68  table. 
These  corrections  are  stored  in  disk  area  REGICO  on  TS44.  This 
capability  has  not  yet  been  transferred  to  the  VAX  11/780,  as  dis¬ 
cussed  below. 

(2)  LOCATION:  This  progam  computes  locations  with  the  addition  of  a 
correlation  matrix. 

(3)  SUCCESS:  Location  method  using  Pn  and  Pg.  Crustal  layer  veloci¬ 
ties  are  modified  to  improve  the  location  estimates  by  successive 

I  determinations. 

(4)  SIMJL:  Location  method  using  Pn  and  Pg.  Crustal  layer  velocities 
are  determined  at  the  same  time  locations  are  computed.  This 
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option  also  allows  locations  using  azimuth  information  only. 

ACLOC  requires  an  extra  disk  area  called  REG1C0.  Since  the  alloca¬ 
tion  of  disk  areas  on  the  VAX  was  not  determined  at  the  time  we  wished 
to  move  these  programs,  the  conversion  of  ACLOC  to  the  VAX  was  delayed 
until  some  later  date.  We  developed  SUCCESS  and  SIMUL  to  locate  events 
with  Pn  and  Pg.  The  option  SIMUL  also  contains  the  option-  to  compute 
locations  with  azimuth  information  obtained  from  the  Smart  processor. 
Therefore,  two  of  the  options,  LOCATION  and  SIMUL,  were  converted  and 
merged.  (The  conversion  of  SUCCESS  was  delayed  for  the  same  reasons  as 
the  conversion  of  ACLOC.) 

2_.2_.  Formation  of  Modules 

Since  all  of  the  subroutines  and  most  parts  of  the  main  routine  in 
each  of  these  different  options  are  identical,  merging  the  options  to 
form  a  single  program  was  highly  desirable.  This  was  not  possible  on 
the  TS44  system,  because  all  the  options  were  stored  as  independent  pro¬ 
grams.  The  convenient  file  structure  on  the  VAX  enabled  such  a  merger 
to  be  easily  accomplished.  First,  the  main  subroutine  LOK  in  both 
options  was  divided  into  several  sections,  so  that  no  duplicate  coding 
was  needed.  Second,  only  those  portions  with  different  coding  were 
filed  under  separate  names.  Third,  a  driver  routine  was  created  to  pick 
up  various  modules  according  to  the  option  desired  for  execution.  This 
method  also  enables  us  to  add  more  options  in  the  future* 

2_.3_.  Directory  of  Station  Coordinates 

The  directory  of  seismic  stations  and  station  coordinates  contain¬ 
ing  approximately  2950  stations  is  called  STATCO,  and  it  is  stored  in 
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the  TS44S0R  disk  area.  Geographic  coordinates  of  stations  in  degrees, 
minutes,  and  tenths  of  seconds  are  stored  in  this  file.  Additional 
information  includes  geographic  and  geocentric  co-latitudes  and  east 
longitudes  of  stations  in  integer  format,  elevation  and  alphanumeric 
station  information.  This  file  has  been  in  use  at  SDAC  for  many  years. 
However,  updating  station  coordinates  has  not  been  done  for  some  time. 
For  example,  station  names  and  coordinates  of  Alaskan  sites,  SRO  and 
ASRO  stations,  etc.,  are  not  stored  in  STATCO. 

The  station  coordinate  file  on  the  VAX  is  called  STATCOORD.DAT. 
This  file  contains  station  co-latitude  and  longitude  in  radians,  and 
some  information  about  the  status,  history,  operating  conditio ns, and 
elevation  of  the  stations.  More  than  3500  stations  are  in  the  file. 

After  making  some  format  changes  to  accommodate  the  new  station 
coordinate  files,  we  have  compared  some  station  coordinates  in  both  sys¬ 
tems.  Some  discrepancies  were  found  in  the  deci-second  field.  However, 
these  discrepancies  were  all  traced  to  errors  in  STATCO  on  the  TS44  sys¬ 
tem. 


£.4.  Arrival  Time  Conversions 

The  arrival  time  inputs  in  the  location  program  are  formatted  in 
integers  of  hour,  minute,  and  deci-seconds.  Upon  converting  them  into 
floating  point  values  in  seconds,  some  discrepancies  at  the  deci-second 
level  were  found  between  the  TS44  version  and  the  VAX  version  of  loca¬ 
tion  programs.  However,  these  errors  do  not  occur  all  the  time.  The 
error  was  caused  by  the  fact  that  the  IBM  computer's  floating  point  for¬ 
mat  does  not  have  enough  bits  to  contain  all  24  hours  of  time  in 
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seconds.  Consequently*  the  least  significant  bit  would  be  lost  on  some 
occasions  in  which  the  time  of  the  day  was  close  to  the  end  of  the  day. 
The  VAX  version  of  time  conversion  has  not  displayed  this  type  of  prob¬ 
lem. 


Matrix  Inversion  Routine 

Matrix  inversion  requires  a  large  memory  area  and  a  long  computing 
time.  It  was  not  possible  to  install  an  accurate  matrix  inversion  rou¬ 
tine  in  the  TS44  system  because  of  its  core  and  computation  time  limita¬ 
tions. 

The  location  method  with  correlation  matrix  was  worked  out  in  a 
previous  paper  (Chang  et  al .  1980).  However*  it  was  not  possible  to 
compute  an  accurate  inverse  matrix  at  that  time.  The  result  of  loca¬ 
tions  with  a  correlation  matrix  was  consequently  tested  with  only  one 
event,  L0NGSH0T.  A  total  of  184  station  inputs  was  divided  into  9  sta¬ 
tion  groups,  and  inversions  were  carried  out  in  the  corresponding  block 
diagonal  form. 

A  robust  matrix  inversion  routine  was  adopted  from  the  IMSL  library 
and  installed  in  the  VAX  location  program.  This  is  a  significant 
achievement  in  that  we  can  now  use  the  location  with  the  correlation 
matrix  option  for  any  event  without  station  groupings.  A  test  run  of 
the  L0NGSH0T  event  with  a  186  by  186  matrix  was  successful.  Since  the 
location  program  is  currently  limited  to  a  maximum  of  200  station 
inputs,  we  think  the  test  shows  that  the  program  can  satisfactorily  han¬ 
dle  the  maximum  input  case.  The  VAX  takes  approximately  2.5  minutes  of 
CPU  time  to  invert  a  matrix  of  this  size. 
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,2 .6,.  Distance  Azimuth  Computation 

Iu  the  past,  source-to-station  distance  and  azinuth  computations 
have  been  carried  out  in  single  precision.  Single  precision  computation 
is  satisfactory  in  most  cases,  but  we  found  that  computations  in  double 
precision  were  necessary  when  two  stations  were  near  each  other.  It  was 
necessary  to  determine  the  effect  of  this  error  on  the  location  program. 

The  subroutine  DIAZ  in  the  ADAPS  system  was  adopted  and  converted 
to  a  subroutine  with  double  precision  computations.  Distance  computa¬ 
tions  between  two  stations  have  been  compared  with  the  existing  subrou¬ 
tine  BJDAZ  and  the  newly  created  double  precision  subroutine  ACDIAZ.  It 
was  found  that  differences  are  very  small  in  most  cases;  all  were  within 
50  meters  of  distances  computed  with  the  double  precision  routine. 
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3.  EVALUATION  OF  TELESEISMIC  LOCATION  METHODS 


2*1*  Ixeai  Locations  ia  Laterally  Bstsiaasasaia  Esrtfr 

Geiger's  least  squares  method  of  locating  seismic  events  determines 
the  event  location  by  minimizing  the  travel  time  residuals,  i.e.,  the 
differences  of  observed  and  theoretical  travel  times  computed  with  a 
standard  Earth  model.  However,  the  Earth  is  in  fact  not  laterally  homo¬ 
geneous.  Therefore,  locations  computed  by  a  least  squares  fit  to  a 
laterally  homogeneous  Earth  model  have  errors  due  to  the  Earth's  hetero¬ 
geneities.  The  amount  of  error  will  depend  on  the  distribution  of  sta¬ 
tions  and  on  the  deviations  of  true  Earth  structures  in  those  station 
areas.  Traditionally,  station  travel  time  corrections  are  added  to  com¬ 
pensate  for  such  errors.  But  it  has  been  found  that  a  simple  constant 
correction  for  each  station  may  not  be  sufficient  for  waves  emerging 
from  all  directions. 

A  simple  method  of  computing  travel  times  for  a  laterally  hetero¬ 
geneous  Earth  model  for  each  source-to-station  path  was  developed  by 
Chang  et  al.  (1980).  In  this  method,  called  ABSLOC,  halfway  travel 
times  of  the  rays  computed  with  Herrin  1968  model  for  various  distance 
ranges  i  and  depth  d,  T(H,d,i),  are  computed  and  subtracted  from  the 
halfway  travel  times  computed  for  other  regional  Esrth  models  A,  B,  C, 
D,  E,  etc.,  thus  resulting  in  tables  of  halfway  travel  time  corrections 
for  those  regional  models  with  respect  to  the  Herrin  model. 

5 t(A,d, i)  -  T(A,d, i)  -  T(H,d,i)  (1) 
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These  halfway  travel  tiwe  corrections  are  stored  in  the  computer  in 
tabular  fora.  One  can  obtain  halfway  travel  tiaes  of  any  aodel  siaply 
by  adding  these  corrections  to  the  halfway  travel  tiae  of  the  Herrin  '68 
aodel,  i.e., 

T(A,d,i)  -  KH.af.i)  ♦  &  t(A,d, i)  (2) 

The  travel  tiae  of  a  laterally  heterogeneous  Earth,  starting  froa  the 
source  in  region  A  and  emerging  at  a  station  in  region  B,  is  equal  to 
the  standard  Herrin  travel  time  for  distance  i  plus  travel  tiae  correc¬ 
tions  for  regional  models  A  and  B, 

T(A,B,d,i)  -  T(H,d,i)  ♦  6t(A,d,i)  +  T(H,0,i)  ♦  <$t(B,0,i) 

-  TH(d,i)  ♦  5 t(A,d, i)  +  6t(B,0,i)  (3) 

where  TH(d,i)  is  the  standard  Herrin  travel  time  at  source  depth  d  and 
distance  i.  Seismic  event  location  with  ray  paths  through  a  laterally 
heterogeneous  Earth  is  made  possible  in  this  manner.  A  test  of  this 
method  using  four  nuclear  explosions  showed  that  errors  were  reduced  to 
about  half  of  the  errors  computed  with  standard  net hods  (loc.  cit.). 
However,  it  was  necessary  to  assign  regional  models  manually  for  each 
source-to-station  path;  furthermore,  only  four  regional  models  were 
used . 

In  the  subsequent  year,  Chang  et  al.  (1981a)  modified  the  method  so 
that  regional  models  for  source  and  station  areas  can  be  assigned 
without  being  input  manually  by  the  analyst.  This  new  method,  called 
ACLOC,  uses  the  729  geographic  regions  of  Flinn  et  al.  (1974),  and  one 
regional  model  is  assigned  to  each  geographic  region.  A  total  of  20 
regional  models  was  adopted  from  the  published  literature  and  computa- 
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tions  were  made  for  travel  tine  correctiona.  However,  alnoat  half  of 
these  models  were  found  to  be  inadequate  because  the  resulting  travel 
time  corrections  did  not  agree  with  published  travel  time  residuals.  By 
eliminating  the  inadequate  models,  a  total  of  11  regional  models  was 
left  in  the  ACLOC  method. 

The  effectiveness  of  the  ACLOC  method  was  compared  against  the 
traditional  location  method  using  Herrin  '68  travel  times.  A  total  of 
12  nuclear  explosions  were  used  for  this  comparison,  and  location  errors 
of  these  methods,  tested  for  depth-free  and  depth-restricted  runs,  are 
shown  in  Table  1.  Results  in  Table  I  show  that  the  ACLOC  method  gives 
better  location  estimates  than  HERRIN,  in  either  depth-free  or  depth- 
restricted  runs.  On  the  average,  the  location  errors  are  reduced  by 
about  3  kilometers. 

One  reason  that  greater  inprovements  in  errors  were  not  observed  is 
that  there  were  many  seismic  stations  that  lie  closer  than  30  degrees 
from  the  event.  Seismic  rays  travelling  through  the  Earth  for  these 
regional  stations  penetrate  only  to  about  800  kilometers,  depending  on 
the  distance  from  the  source.  Since  the  Earth's  crust  and  upper  mantle 
are  known  to  be  quite  heterogeneous,  these  stations  are  the  ones  that 
are  most  influenced  by  incorrect  regional  models.  Subsequent  tests, 
using  only  stations  at  distances  greater  than  24  degrees,  resulted  in 
slightly  better  locations.  These  results  are  shown  in  Table  II. 

The  ACLOC  method  should  give  better  results  because  it  corrects  for 
lateral  heterogeneities.  The  method  is  limited  by  our  lack  of  knowledge 
of  the  Earth's  heterogeneity,  which  appears  in  two  ways.  The  first  is 
that  the  amdels  we  choose  are  not  adequate;  and  secondly,  the  selection 
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of  the  model  assigned  to  some  of  the  areas  may  not  be  adequate.  The 
adequacy  of  the  assigned  models  was  examined  by  comparing  travel  times 
computed  through  them  to  published  travel  time  residuals.  Further ,  one 
can  calibrate  some  of  the  model  assignments  by  comparing  the  published 
residuals  of  some  stations  against  the  corrections  assigned  by  the  com¬ 
puter.  This  was  done  by  checking  the  residuals  of  the  ACLQC  runs  and 
re-assigning  regional  models  to  those  stations  that  had  high  residuals. 
The  calibrated  version  is  called  ACL0C3. 

In  addition  to  the  travel  time  residuals  data  published  by  Sengupta 
and  Julian  (1974)  and  by  Cleary  and  Hales  (1966),  we  acquired  a  copy  of 
the  worldwide  travel  time  residuals  published  by  Poupinet  (1979).  This 
latter  paper  covers  more  worldwide  areas  than  the  previous  papers. 
Thus,  another  version  called  ACIXJC4,  with  regional  model  assignments 
based  on  Poupinet's  residuals,  was  created.  Table  ITT  compares  the 
location  errors  of  these  various  versions,  HERRIN,  ACLQC,  ACL0C3,  and 
ACL0C4.  The  results  of  ACL0C4  are  about  the  same  as  those  of  ACLOC;  and 
the  results  of  ACLOC 3,  which  used  the  residuals  of  Sengupta  and  Julian 
(1974)  and  of  Cleary  and  Hales  (1966),  are  the  best  among  these  ver¬ 
sions.  It  is,  however,  possible  that  ACL0C4  will  give  better  locations 
when  the  model  assignments  are  corrected  for  bad  stations.  The  average 
improvement  in  location  errors  as  compared  to  HERRIN  in  Table  III  is 
about  6  kilometers , declining  from  21  to  15  kilometers. 
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TABLE  I 


Comparison  of  Location  Errors.  ACLOC  vs.  HERRIN 
(Stations  within  12  degrees  are  excluded) 


Event 

Name 

No. 

Sta. 

HERRIN-R 

Location  Errors 
ACLOC-R 

in  Kilometers 
HERRIN-F 

ACLOC-F 

L0NGS1UT 

184 

11.04 

7.70 

19.42 

16.60 

GNOME 

44 

43.83 

45.68 

43.83 

46.75 

SHOAL 

21 

51.64 

24.88 

41 .52 

44.07 

SALMON 

75 

22.91 

16.36 

20.73 

16.15 

CORDUROY 

89 

12.19 

2C  08 

20.93 

20.80 

GASBUGGY 

60 

21.19 

9.00 

31.30 

22.18 

FAULTLESS 

106 

6.98 

7 .46 

5.46 

2.46 

PURSE 

121 

9.77 

9.27 

14.89 

12.45 

RULXSON 

74 

12.88 

/  .08 

15.51 

10.64 

PIPKIN 

103 

10.32 

7.76 

26.05 

22.96 

TERRI NE 

53 

6.03 

12.44 

24.05 

16.44 

CARPETBAG 

136 

7.C9 

8.54 

11.66 

7.41 

Average 

17.99 

14.69 

22.95 

19.91 

R  -  depth  restricted,  F  -  depth  free 
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TABLE  II 


Comparison  of  Location  Errors(  ACLOC  vs.  HERRIN 
(Stations  within  24  degrees  are  excluded) 


Event 

Name 

No. 

Sta. 

HERRIN-R 

Location  Errors 
ACLOC-R 

in  Kilometers 
HERRIN-F 

ACLOC-F 

LONG SHOT 

168 

21.85 

18.05 

21.32 

17.10 

GNOME 

11 

15.93 

34.47 

56.31 

59.48 

SHOAL 

10 

16.66 

16.56 

42.29 

44.81 

SALMON 

26 

17.64 

8.49 

14.68 

7.58 

CORDUROY 

72 

9.64 

6.56 

8.88 

7.39 

GASBUGGY 

40 

22.43 

21.74 

31.77 

27.63 

FAULTLESS 

101 

5.58 

2.01 

4.66 

2.19 

PURSE 

108 

13.98 

11.76 

16.88 

13.07 

RULISON 

55 

1.86 

4.09 

6.51 

4.19 

PIPKIN 

85 

16.57 

15.41 

20.21 

15.26 

TERRI NE 

42 

16.92 

15.81 

19.08 

13.66 

CARPETBAG 

118 

15.32 

13.12 

15.80 

13.30 

Average 

14.53 

14.00 

21.53 

18.81 

R  -  depth  restricted,  F  -  depth  free 
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A  comparison  of  the  depth  errors  for  various  methods  is  also  shown 
in  Table  III.  On  the  average  we  find  that  depth  errors  are  two  to  three 
tines  larger  than  location  errors;  this  is  a  well  known  phenomenon  for 
shallow  events.  If  we  compare  the  average  depth  errors  of  HERRIN  to 
those  of  ACLOC  or  ACL0C3,  we  see  that  depth  errors  also  are  reduced  in 
the  ACLOC  method.  However,  without  accurate  knowledge  of  regional  geo¬ 
logic  structures,  further  improvements  of  the  ACLOC  method  are  rather 
difficult. 

It  is  noteworthy  that  all  depths  computed  with  the  HERRIN  method 
are  positive,  while  ACLOC  and  ACL0C3  methods  resulted  in  some  negative 
depths.  We  also  note  that  those  events  with  negative  depths  are  all  at 
NTS.  Since  the  mean  depth  over  an  ensemble  of  events  at  NTS  should  be 
near  zero,  or  perhaps  1  km,  this  suggests  that  the  source  model  assign¬ 
ments  for  NTS  may  be  proper. 
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TABLE  III 


Comparison  of  Three  Versions  of  ACLOC  vs.  HERRIN 
(Stations  within  24  degrees  are  excluded,  all  runs  depth  free) 


Event  Location  and  Depth  Errors  in  Kilometers 


Marne 

HERRIN  (d) 

ACLOC  (d) 

ACLOC 3  (d) 

ACL0C4  (d) 

LONGS HOT 

21.32  (98.9) 

17.10  (88.0) 

13.53  (92.0) 

14.94  (77.2) 

GNOME 

56.31  (89.4) 

59.48  (37.8) 

39.41  (37.8) 

58.30  (80.6) 

SHOAL 

42.29(165.2) 

44.81(158.8) 

30.83(158.8) 

46.92(163.7) 

SALMON 

14.68  (35.4) 

7.58  (12.6) 

12.95  (12.6) 

7.03  (35.9) 

CORDUROY 

8.88  (27.8) 

7.39(-13.3) 

10.57 (-18.9) 

9.96  (35.8) 

GASBUGGY 

31.77  (88.2) 

27  .63  (  55.8) 

19.71  (50.2) 

24.15  (84.5) 

FAULTLESS 

4.66  (33.7) 

2.19(-24.4) 

3.50  (53.7) 

6.71  (30.7) 

PURSE 

16.88  (57.5) 

13.07  (28.2) 

7.77  (28.7) 

13.98  (54.3) 

RULISON 

6.51  (73.9) 

4.19  (38.3) 

4.27  (41.9) 

8.42  (71.6) 

PIPKIN 

21.21  (36.6) 

15.26  (-1.7) 

17.28  (-2.4) 

20.57  (37.2) 

TERRI NE 

19.08  (21.8) 

13.66(-23.2) 

11.16(-24.3) 

12.00  (13.3) 

CARPETBAG 

15.80  (51.1) 

13.30  (23.3) 

6.92  (23.0) 

12.77  (50.0) 

Average 

21.53  (65.0) 

18.81  (42.1) 

14.83  (43.4) 

19.65  (61.2) 

*  Average  depth  errors  are  the  average  of  absolute  values. 
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3. .2_ .  Locations  With  Correlation  Matrix 

The  ACLOC  method  discussed  in  the  previous  section  assumes  that 
travel  time  residuals  are  fixed  deterministic  quantities  that  can  be 
independently  estimated  for  each  source  and  station.  Alternately,  one 
can  eliminate  the  location  bias  due  to  clusters  of  seismic  station 
located  in  a  close  proximity  by  weighting. 

The  second  method  of  locating  teleseismic  events,  LOCATION,  assumes 
that  stations  in  close  proximity  may  be  measuring  the  same  effect  of 
regional  structure.  Thus  by  inserting  a  correlation  matrix  in  the  nor¬ 
mal  equation,  the  effect  of  clustered  stations  will  be  corrected  by  the 
weighting  of  the  travel  time  residuals  in  accordance  with  their  overall 
intercorrelation  matrix.  The  advantage  of  this  method  is  that  the 
correction  of  station  bias  improves  not  only  the  location  estimates,  but 
also  the  confidence  intervals,  which  will  not  contract  substantial ly  as 
arrival  times  from  highly  correlated  stations  are  added  to  the  data 
base. 

In  the  previous  work,  correlation  coefficients  of  residuals  were 
estimated  by  using  186  station  residuals  of  LONGSHOT  event  (Chang  et  al. 
1980).  Coefficients  of  0.7  for  0.1°<  £<  5°,  and  0.3  for  5°  <  A  <  10° 

were  used  to  test  the  method.  However,  the  method  was  not  thoroughly 
tested  because  an  accurate  matrix  inversion  routine  was  not  available  on 
the  TS44  system.  The  new  version  of  LOCATION  on  VAX  is  capable  of  fully 
testing  this  method. 

The  correlation  coefficients  on  the  VAX  location  program  are  com¬ 
puted  with  the  following  equations: 


c  =  0.7  exp(-( A/7 °) 2)  ;  for  A  >  0  degree,  and 


c  =  1 .0  ;  for  A  *  0  degree. 

Values  computed  with  this  equation  are  in  agreement  with  those  used  in 
the  previous  work,  but  are  better  in  the  sense  that  the  values  are  con¬ 
tinuous  with  distance.  However,  correlation  structures  should  be 
further  investigated  in  the  future.  The  question  of  whether  these 
correlation  coefficients  vary  from  one  region  to  another  is  also  an 
important  matter  to  be  determined  in  the  future. 

A  total  of  15  nuclear  explosions  were  used  to  evaluate  the  U)CAT10N 
method  of  locating  events  with  correlation  matrix.  These  15  events 
include  12  explosions  used  in  Tables  I  -  III,  plus  two  French  Sahara 
tests  and  a  Kazakh  event.  The  location  coordinates  of  the  French  tests 
were  given  by  Bolt  (1976).  The  Kazakh  event  was  a  cratering  event,  and 
the  location  of  the  crater  was  determined  from  a  satellite  photograph 
(Fitch  and  North,  1980). 

Comparison  of  location  errors  using  a  standard  location  method  with 
the  Herrin  table  are  compared  with  errors  using  the  correlation  matrix. 
The  result  shows  that  better  locations  can  be  obtained  using  the  corre¬ 
lation  matrix  method.  However,  improvements  are  not  as  great  as  for  the 
LONGSHDT  event  in  most  cases.  The  reason  for  this  is  that  for  some 
events  like  LONGSHOT,  many  stations  were  situated  in  close  proximity; 
thus,  for  these  events,  station  bias  is  removed  by  using  the  correlation 
matrix  method.  For  the  other  events,  stations  are  not  distributed  in 
close  proximity  to  each  other,  thus  stations  exerting  regional  biases 
are  not  overweighted ,  and  there  are  no  appreciable  improvements.  One 
can  compare  these  events  in  Table  IV  with  those  corresponding  events 


-23- 


with  the  ACLOC  method  in  Table  I,  and  find  that  the  the  improvement  is 
about  the  same.  Note  that  both  methods  are  better  than  standard  loca¬ 
tion  results. 

The  Kazakh  event  is  the  only  event  for  which  SRST  corrections  are 
available  (Veith,  1974).  We  have  relocated  this  event  in  two  ways:  with 
the  SRST  corrections,  and  without  them.  Locations  were  sightly  better 
when  the  SRST  corrections  were  applied  to  the  HERRIN  location  method, 
but  no  improvements  were  observed  when  SRST  corrections  were  used 
together  with  the  correlation  matrix  method.  In  the  latter  case,  the 
lack  of  improvement  using  SRSTs  may  be  attributable  to  the  fact  that  a 
number  of  stations  used  for  the  Kazakh  event  were  not  included  in  the 
SRST  corrections  table. 

Absolute  depth  errors  for  the  1OCATI0N  method  are  about  the  same  as 
depth  errors  for  the  ACLOC  method.  The  advantage  of  the  LOCATION  method 
over  the  ACLOC  method  is  that,  without  knowing  the  details  of  the  bias¬ 
ing  geology,  it  minimizes  the  effects  of  differences  between  ray  path* 
in  the  upper  mantle.  This  al lows  the  LOCATION  method  to  use  all  the 
station  inputs.  However,  if  we  compare  the  results  of  ACLOC  methods 
using  all  stations  (Table  I),  we  find  that  the  difference  in  location 
errors  are  very  small. 

If  the  station  bias  is  properly  removed  by  the  LOCATION  method,  the 
shape  of  the  resulting  error  ellipse  should  also  be  more  symmetrical. 
Removal  of  bias  will  also  increase  the  chance  of  containing  the  true 
location  within  the  computed  error  ellipse.  Test  results  in  Table  IV 
confirm  this.  In  spite  of  only  slight  improvement  in  location  errors, 
more  events  relocated  by  the  LOCATION  method  resulted  in  error  ellipses 


including  the  true  location.  In  HERRIN-F  runs,  only  two  out  of  16 
events  contain  the  true  locations  within  the  computed  error  ellipses. 
For  two  of  the  remaining  events,  the  true  locations  lay  outside  of  the 
three-dimensional  ellipsoids,  but  within  the  two-dimensional  ellipses; 
and  for  twelve  events,  the  true  locations  lay  outside  of  both  the  three- 
and  two-dimensional  error  ellipses.  In  contrast,  there  were  only  four 
out  of  16  events  in  LOCATION-F  runs  for  which  two-dimensional  error 
ellipses  did  not  contain  the  true  locations,  and  five  events  where  the 
true  locations  lie  outside  both  the  three-  and  two-dimensional  ellipses. 
One  event  (RULISON)  resulted  in  a  true  location  inside  the  three- 
dimensional  ellipsoid,  but  outside  the  two-dimensional  error  ellipse. 
This  unusual  result  will  be  discussed  in  the  appendix. 

Better  epicenters  were  obtained  when  event  depths  were  restricted, 
with  only  one  event  in  HERRIH-R  runs  not  ntaining  the  true  location  in 
the  error  ellipse.  For  three  out  of  six  LOCATION-R  runs,  the  true  loca¬ 
tions  were  inside  the  error  ellipses. 
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TABLE  IV 


Comparison  of  Location  Errors*  LOCATION  vs.  HERRIN 
(Stations  within  12  degrees  are  excluded) 


Event  Location  and  Depth  Errors  in  Kilometers 


Name 

HERRIN-R 

LOCATION- R 

HERRIN-F 

(d) 

LOCATION- F 

(d) 

LONGSHOT 

611.29 

~T 

4.90  | 

#19.51 

1 

(82.0) | 

*10.76 

(65.2) 

GNOME 

643.78 

642.00  1 

#43.78 

(-0.6) 1 

#42.07 

(3.1) 

SHOAL 

51.73 

37.01  1 

*41 .44 

(115.4)  1 

41.17 

(109.1) 

SALMON 

622.81 

12.97  1 

#20.68 

(20.3) 1 

12.42 

(-23.7) 

CORDUROY 

12.31 

1.57  1 

#21.03 

(72.6)  I 

#23.08 

(80.0) 

GASBUGGY 

621 .46 

630.68  | 

#31.51 

(40.5) | 

#35.12 

(32.7) 

FAULTLESS 

7.03 

5.54  | 

*  5.46 

(49.3) | 

*5.92 

(55.1) 

PURSE 

6  9.94 

7.37  | 

#14.97 

(56.5)  I 

#15.19 

(58.7) 

RULISON 

613.11 

614.49  1 

#15.66 

(22.4) I 

$15.78 

(13.2) 

PIPKIN 

10.38 

16.83  1 

#26.05 

(67.6) 1 

#30.26 

(69.6) 

TERRI NE 

5.80 

4.71  I 

#24.21 

(64.7) 1 

*23.52 

(68.9) 

CARPETBAG 

7.18 

6.87  1 

#11.73 

(61.6)  1 

*12.80 

(57.3) 

SAPHIR 

2.00 

2.39  1 

2.54 

(-4.1) 1 

2.53 

(3.3) 

RUBIS 

9.11 

6.61  1 

8.49 

(-5.8)  1 

5.94 

(-6.4) 

KAZAKH 

610.55 

5.60  1 

#10.49 

(0.9)  1 

5.45 

(2.9) 

KAZAKH(SRST) 

6  6.79 

5.71  1 

#  8.13 

(-17.0) { 

6.52 

(-15.3) 

Average 

15.33 

1 

12.83  1 

19.11 

(42.6) 1 

18.03 

(41.5) 

R  -  Depth  restricted 
F  -  Depth  free 

*  -  True  location  lies  outside  of  the  3D  95Z  F-statistic  ellipse 

#  -  True  location  lies  outside  of  the  3D  and  2D  95Z  F-statistic  ellipses 

&  -  True  location  lies  outside  of  the  2D  95Z  F-statistic  ellipse 

$  -  True  location  lies  inside  of  the  3D*  but  outside  of  the  2D 

95Z  F-statistic  ellipses. 


-26- 


4.  SEISMIC  EVENT  LOCATIONS  WITH  REGIONAL  DATA 


4.1.  Locations  with  Pn .  Pg  ind  Lg. 

The  standard  method  of  locating  events,  at  both  regional  and 
teleseismic  distances,  computes  event  locations  with  travel  times  com¬ 
puted  from  a  standard  earth  model*  Because  the  layer  velocities  of  the 
crust  and  upper  mantle  can  change  in  a  rather  short  distance,  location 
accuracies  using  Pn  and  Pg  may  vary  from  region  to  region.  It  is  known 
that  better  locations  may  be  obtained  if  an  adequate  regional  model  is 
used,  but  such  a  regional  model  is  not  available  for  most  areas. 

Methods  developed  in  the  previous  year  are  designed  to  compensate 
for  such  local  variations  by  modifying  layer  velocities  while  computing 
event  locations.  Details  of  these  methods  are  discussed  in  Chang  et  al. 
(1981b). 

On  the  VAX  version  of  the  simxltaneouB  inversions  method  (SIMUL), 
the  capability  of  using  Lg  arrivals  was  added  to  the  capability  of  using 
Pn  and  Pg.  The  Lg  phase  is  assumed  to  have  an  initial  velocity  of  3.54 
km/sec,  and  this  value  is  modified  during  the  iteration  cycle.  It  is 
normally  difficult  to  use  Lg  for  locating  events  because  its  velocity  is 
not  well  determined.  However,  the  SIKJL  method  can  overcome  this  diffi¬ 
culty  by  iterative  modifications. 

We  will  now  compare  the  location  errors  using  Pn,  Pg  (SIMUL)  with 
the  errors  from  using  Pn,  Pg  and  Lg  (PLUSLG)  in  application  to  data  from 
12  explosions. 
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The  location  program  on  the  VAX  was  programmed  with  Julian's  (1973) 
original  idea  of  weighing  each  arrival  datum.  In  order  to  use  all  the 
regional  data  for  each  event  in  a  single  location,  errors  of  observation 
of  1.0  second  for  Pn  and  3.0  seconds  for  Pg  and  Lg  are  adopted  from 
McCowan  and  Needham  (1978).  The  intended  effect  of  this  was  to  reduce 
the  dependence  of  the  epicenter  on  those  arrival  times  which  cannot  be 
read  as  accurately  as  Pn.  Each  observation  (diagonal  element  of  the 
covariance  matrix)  is  multiplied  by  the  inverse  of  these  error  values, 
so  that  Pg  and  Lg  are  weighted  as  only  one  third  of  the  Pn  observation. 
However,  we  are  not  certain  that  these  error  values  are  suitable  for 
general  use*  In  order  to  show  the  effect  of  weighting,  the  12  explo¬ 
sions  were  re-located  in  two  ways;  one  as  above,  and  the  other  weighting 
all  phases  equally.  Comparison  of  the  location  errors  are  shown  in 
Table  V. 

Many  interesting  conclusions  can  be  drawn  by  examining  the  results 
in  Table  V.  These  are: 

(1)  The  addition  of  Lg  arrivals  in  locating  events  did  not  improve  the 
locations.  However,  no  significant  deterioration  is  seen  either* 
The  usefulness  of  Lg  phase  in  location  projects  depends  on  how 
accurately  these  arrival  times  can  be  determined,  and  whether  there 
are  azimuthal  variations  in  Lg  velocities. 

(2)  Changing  the  weighting  factors  did  not  change  the  location  errors 
appreciably.  However,  changing  these  weighting  factors  to  seme 
other  values,  for  example,  1.5  for  Pg  and  3.0  for  Lg,  may  be  tried 
in  the  future. 
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(3)  Table  V  shows  12  events  which  were  relocated  four  different  ways. 
For  each  of  the  four  experiments  location  errors  are  mostly  within 
10  kilometers  of  the  true  location.  In  each  experiment  the  average 
error  is  about  6  kilometers.  This  result  is  significantly  better 
than  location  errors  with  using  teleseismic  data,  as  shown  in 
Tables  I  through  IV. 

(4)  Although  there  are  a  number  of  input  data  in  each  event,  actually 
there  were  only  a  few  stations  used  in  these  events.  By  examining 
the  number  of  Pn.  Pg  and  Lg  used  in  each  case,  one  finds  that  the 
number  of  stations  in  these  events  ranges  from  4  to  13.  Although 
this  number  may  be  typical  for  most  seismic  events,  we  point  out 
that  no  appreciable  improvements  in  Pn.  Pg  and  Lg  velocities  can  be 
made  with  only  four  or  five  stations.  With  this  small  number  of 
stations  it  is  remarkable  that  location  errors  with  regional  data 
are  much  better  than  teleseismic  locations. 

(5)  With  the  exception  of  GASBUGGT,  better  locations  can  be  obtained 
with  an  increased  number  of  stations. 

(6)  Examination  of  residuals  in  GASBUGG?  shows  that  4  out  of  12  Pg 
inputs  show  large  residuals,  and  2  out  of  6  Lg  residuals  are  also 
high.  TVo  bad  cases  of  Lg  residuals  suggest  that  Lg  velocities  may 
be  different  in  different  directions  of  propagation.  However,  it 
is  interesting  to  note  that  although  bad  residuals  are  observed  in 
Pg  and  Lg  arrivals,  location  errors  are  better  when  all  phases  are 
weighted  equally. 


TABLE  V 


Comparison  of  Location  Errors,  SIM7L  vs.  PLUSLG 
(Depth  restricted) 


Event 

Marne 

No. 

Pn 

No. 

Pg 

No. 

I* 

Unequal  Weights 
SIMJL  PLUSLG 

Equal  Weights 
SIMJL  PLUSLG 

FAULTLESS 

5 

5 

5 

*18.65 

11.99 

7.37 

1.19 

RULISON 

10 

7 

4 

0.69 

0.75 

5.11 

3.48 

PASSAIC 

4 

3 

4 

3.04 

3.07 

3.49 

*  3.89 

ROCKVILLE  DAM 

9 

8 

8 

4.43 

3.78 

16.02 

*10.69 

DORMOUSE  PRIME 

3 

4 

4 

7.00 

7.84 

10.21 

8.58 

KLICKITAT 

13 

12 

12 

0.76 

1.18 

0.35 

1.49 

BANDICOOT 

5 

4 

5 

8.66 

2.23 

2.23 

6.50 

SHOAL 

12 

0 

4 

1.88 

1.49 

1.88 

4.31 

MERRIMAC 

9 

0 

8 

5.77 

7.11 

5.77 

14.57 

GAS BUGGY 

12 

12 

6 

*10.03 

*  9.68 

*  8.92 

*  8.54 

PILEDRIVER 

5 

5 

2 

5.73 

5.73 

4.53 

4.53 

ROANOKE 

4 

4 

4 

4.10 

7.05 

5.08 

*10.28 

Average 

5.90 

5.16 

5.91 

6.50 

*  -  True  location  lies  outside  of  the  2D  95X  E-statistic  ellipse 
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4.2,.  Location  with  All  Phases 

So  far  we  have  been  evaluating  location  methods  with  either  telese- 
ismic  P  phases  or  with  regional  phases.  These  two  methods  were 
separated  during  the  development  stage,  because  they  present  different 
problems  in  locating  seismic  events.  However,  since  events  are  located 
with  all  available  phases,  the  different  versions  must  be  combined  and 
the  techniques  evaluated  with  all  phases. 

Two  methods,  LOCATION  and  SIMUL,  were  put  together  so  that  this 
program  will  compute  locations  with  all  phases.  In  this  final  program 
all  regional  phases  use  the  similtaneous  inversions  method  to  modify  Pn, 
Pg  and  Lg  velocities  during  the  iterative  least  squares  cycles,  and  a 
matrix  is  inserted  in  the  normal  equation  so  that  the  correlation  matrix 
method  can  be  performed.  If  one  chooses  not  to  use  the  correlation 
matrix  method,  a  unit  matrix  consisting  of  one's  on  the  diagonal  ele¬ 
ments  and  xero's  elsewhere  is  substituted  for  the  correlation  matrix,  so 
there  is  no  change  in  computing  the  solutions.  Otherwise,  non-zero 
coefficients  will  be  assigned  to  each  elements  as  described  in  section 
3.2. 

The  non-diagonal  elements  of  regional  phases  are  all  assigned  zero 
values.  This  was  necessary  because  we  do  not  know  proper  values  for 
them.  However,  correlation  coefficients  for  regional  phases  should  be 
determined  in  the  future. 

Six  explosions  were  used  to  compare  location  errors  without 
(NOMATX)  and  with  a  correlation  matrix  (MATRIX).  Results  are  shown  in 
Table  VI.  Location  errors  using  either  method  are  as  high  as  those  with 
teleseismic  P  only,  see  Table  IV.  Improvement  of  errors  using  the 
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MATRIX  method  as  compared  to  the  NOMATX  method  is  about  4  kilometers 
from  26  to  22  kilometers.  Using  the  correlation  matrix  method,  only  one 
event,  GASBUGGY,  resulted  in  inferior  locations  as  compared  to  NOMATX. 

Location  errors  given  in  Table  VI  average  to  greater  than  20  kilom¬ 
eters.  This  was  probably  due  to  the  high  weighting  factor  imposed  on 

teleseismic  P.  In  Julian's  (1973)  location  program,  each  station  is 

weighted  by  the  inverse  of  the  assigned  standard  error.  In  the  experi¬ 
ments  shown  in  Table  VI,  standard  errora  are  0.2  for  P,  1.0  for  Pn  and 

3.0  for  Pg  and  Lg.  Consequently,  P  signals  were  weighted  five  times  as 

large  as  those  of  Pn  and  15  times  as  large  as  Pg  or  Lg.  We  think  this 
weighting  is  not  satisfactory.  Although  the  values  of  standard  errors 
for  regional  phases  were  adopted  from  McCowan  and  Needham  (1978),  the 
value  for  P  should  be  about  equal  to  the  value  of  Pn.  Results  of 
weighted  and  unweighted  runs  shown  in  Table  V  indicate  that  locations 
with  equal  weighting  for  all  regional  phases  are  only  slightly  worse. 
This  suggests  that  these  weighting  ratios  needs  to  be  further  investi¬ 
gated.  Note  that  using  a  ■  0.2  second  for  P  leads  to  x2-error  ellipses 
that  do  not  include  the  true  epicenters,  so  cj  should  be  larger  than  0.2 
second  for  P. 

In  Table  VII  we  show  location  errors  of  NOMATX  and  MATRIX,  with 
equal  weighting  factors  assigned  to  all  phases.  The  improvement  in 
locations  are  obvious.  Furthermore,  note  that  errors  in  MATRIX  runs  are 
much  better  than  errors  in  NOMATRIX  runs.  From  this  result  we  conclude 
that  equal  weighting  factors  or  more  weighting  for  regional  phases  may 
be  better.  However,  we  emphasize  that  the  mutual  relations  of  weighting 
factors  for  various  phases  in  the  location  program  should  be  further 
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investigated 


Depth  estimates  of  these  events  are  also  shown  in  Table  VI  and  VII. 
Depths  computed  with  all  phases  are  not  much  different  from  those  depths 
computed  with  only  teleseismic  P  (Table  IV).  This  is  expected  because 
the  simultaneous  inversions  method  does  not  estimate  event  depths.  The 
method  computes  the  best  fit  linear  slopes  for  Pn,  Pg  and  Lg;  but  does 
not  attempt  to  estimate  the  event  depth.  Consequently  the  computed 
depth  using  all  phases  will  resemble  the  depth  computed  with  P  only. 
Improvements  in  depth  estimates  should  be  a  research  topic  for  the 
future. 

True  locations  of  these  events  relocated  with  the  MATRIX  method 
were  all  found  to  lie  within  the  computed  error  ellipses.  On  the  other 
hand,  true  locations  of  events  relocated  with  the  NOMATX  method  were 
mostly  found  outside  of  the  computed  error  ellipses.  When  all  phases 
were  weighted  equally  in  the  NOMATX  method,  more  true  locations  were 
found  within  the  ellipses.  This  was  the  result  of  the  better  location 
obtained  with  the  help  of  regional  phases. 
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TABLE  VI 


Comparison  of  Location  Errors.  NOMATX  vs.  MATRIX 
Using  All  Phases  and  Strong  Weights  on  P 


Event  No.  No.  No.  No.  Depth  Restricted  Depth  Free 


Name 

P 

Pn 

Pg 

Lg 

1 

1 

NOMATX 

MATRIX 

NOMATX 

(depth) 

MATRIX 

(depth) 

FAULTLESS 

106 

5 

5 

5 

T 

1 

1 

6.89 

5.77 

*  5.36 
(50.3) 

*  5.73 
(55.4) 

RULISOH 

74 

10 

7 

4 

1 

1 

612.36 

619.15 

#14.97 

(27.3) 

#19.87 

(11.2) 

SHOAL 

21 

12 

0 

4 

1 

1 

644.32 

28.02 

*25.11 

(102.8) 

*18.05 

(92.1) 

GAS BUGGY 

60 

12 

12 

6 

1 

1 

626.24 

637.33 

#33.68 

(39.2) 

#38.52 

(16.7) 

SALMON 

76 

4 

0 

0 

1 

1 

622.71 

12.86 

#20.54 

(-21.9) 

12.39 

(-23.8) 

GNOME 

44 

18 

1 

0 

1 

1 

_1_ 

641.92 

630.02 

#41.92 

(-0.6) 

#30.02 

(-12.2) 

Average 

1 

1 

25.74 

22.19 

23 .60 

20.76 

*  -  True  location  lies  outside  of  the  3D  95t  F-statistic  ellipse 


#  -  True  location  lies  outside  of  the  3D  and  2D  952  F-statistic  ellipses 
&  -  True  location  lies  outside  of  the  2D  952  F-statistic  ellipse 
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TABLE  VII 


Comparison  of  Location  Errors.  NOMATX  vs.  MATRIX 
Using  All  Phases  and  Equal  Weights  to  All  Phases 


Event  Depth  Restricted  Depth  Free 


Name 

NOMATX 

MATRIX 

NOMATX 

(d) 

MATRIX 

(d) 

FAULTLESS 

3.37 

1.47 

*  2.38 

(47.4) 

*  0.78 

(55.6) 

RULISON 

2.09 

3.41 

2.42 

(17.9) 

3.58 

(5.8) 

SHOAL 

9.83 

1.57 

*  3.92 

(86.7) 

*  3.34 

(82.4) 

GASBUGGY 

&10.11 

&  9.86 

#  9.59 

(29.5) 

S  9.92 

(3.7) 

SALMON 

&21.36 

13.02 

#19.13 

(-22.0) 

12.26 

(-24.6) 

GNOME 

&25.50 

10.77 

#26.95 

(-18.5) 

11.36 

(-4.0) 

Average 

12.04 

6.68 

10.73 

(37.0) 

6.87 

(29.4) 

*  -  True  location  lies  outside  of  the  3D  95%  F-statistic  ellipse 

#  -  True  location  lies  outside  of  the  3D  and  2D  95X  F-statistic  ellipses 

&  -  True  location  lies  outside  of  the  2D  95X  F-statistic  ellipse 

$  -  True  location  lies  inside  of  the  3D,  but  outside  of  the  2D 

95%  F-statistic  ellipses. 


-35- 


5.  CONCLUSIONS  AND  RECOMMENDATIONS 


Seismic  locations  using  teleseismic  P  data  and  using  regional  data 
present  different  problems.  For  this  reason,  our  research  in  the  previ¬ 
ous  years  had  different  location  programs  for  different  phases.  As  a 
result  we  have  developed  many  useful  methods  for  different  problems. 
The  RELS  location  program  on  the  VAX  11/780  at  SDAC  is  a  unified  loca¬ 
tion  program  that  can  handle  regional  and/or  teleseismic  and/or  back 
azimuth  data.  This  program  can  be  executed  with  or  without  the  sinul- 
taneous  inversions  method,  and  with  or  without  the  correlation  matrix 
method.  Since  all  of  these  methods  are  independent  programs  on  the  TS44 
system,  the  unified  location  program  on  the  VAX  is  far  more  versatile 
than  those  separate  programs  on  the  TS44.  A  remaining  task  is  to  add  the 
capability  of  using  source  and  receiver  models  when  the  disk  file  struc¬ 
ture  for  the  VAX  is  agreed  upon. 

Three  location  methods,  two  using  teleseismic  data  and  one  using 
regional  phases,  are  evaluated  by  comparing  location  errors  to  the 
corresponding  errors  computed  by  the  standard  method.  All  methods 
showed  some  improvement  over  the  standard  method;  and  when  all  signals 
were  combined,  resultant  location  errors  were  about  6  kilometers  on  the 
average.  Average  errors  with  teleseismic  P  data  alone  are  much  larger, 
the  best  result  being  for  the  ACL0C3  method,  with  an  average  of  about  14 
kilometers.  Locations  with  regional  data  gave  better  results,  6  kilome¬ 
ters  on  the  average.  These  errors  are  approximately  equal  to  those 
using  al  1  data. 
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Since  the  methods  evaluated  here  do  not  use  travel  time  correc¬ 
tions,  these  methods  are  not  dependent  on  specially  calibrated  stations. 
These  methods  also  do  not  require  a  priori  knowledge  of  the  local  struc¬ 
tures,  therefore  they  can  be  used  to  locate  seismic  events  in  any  parts 
of  the  world.  Locations  computed  with  these  methods  are  absolute  loca¬ 
tions,  not  relative  locations  as  are  those  using  travel  time  corrections 
derived  from  a  master  event. 

The  correlation  matrix  method  was  not  applied  to  the  regional  data, 
because  correlation  coefficients  for  Pn.  Pg  and  Lg,  obtained  at  proxi¬ 
mate  stations  and  the  correlation  coefficients  for  Pn  and  Pg ,  Pn  and  Lg, 
etc.,  observed  at  the  same  station  are  not  known  at  the  present.  Deter¬ 
mination  of  these  coefficients,  and  also  determination  of  standard 
errors  to  be  used  for  weighting  each  input  datum,  should  be  investigated 
in  the  near  future. 

The  RELS  location  program  is  also  capable  of  using  Lg  arrivals  for 
locating  seismic  events.  In  this  study,  all  Lg  arrivals  were  picked  by 
analysts  at  SDAC  at  the  beginning  of  the  Lg  wave  train,  rather  than  at 
the  place  where  the  maximum  amplitude  occurs.  We  also  note  that 
research  on  Lg  has  been  concentrated  on  the  attenuation  phenomena,  but 
not  on  the  regional  variations  of  Lg  velocities.  With  the  additional 
capability  of  locating  events  with  Lg,  the  propagation  phenomena  of  Lg 
in  relation  to  the  earth's  heterogeneities  should  be  further  investi¬ 
gated. 

Although  the  teleseismic  locations  examined  here  use  a  large  number 
of  observations,  the  accuracy  of  location  with  the  teleseismic  data  i6 
worse  than  the  accuracy  using  a  smaller  number  of  regional  data.  There 
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may  be  several  explanations  for  this.  The  first  reason  is  that  dT/dA 
for  regional  phases  are  larger  than  those  of  teleseismic  P,  thus  the 
least-squares  method  of  fitting  a  standard  travel  time  curve  works 
better  for  regional  phases. 

Even  so  we  are  still  surprised  that  one  can  obtain  better  location 
estimates  by  using  data  from  a  handful  of  regional  stations  rather  than 
using  a  few  dozens  of  teleseismic  observations. 

It  is  generally  accepted  that  the  Earth  is  more  heterogeneous  in 
the  crust  and  upper  mantle.  If  so,  the  effect  of  lateral  heterogeneity 
would  appear  stronger  on  the  travel  time  residuals  of  Pn  and  Pg  as  com¬ 
pared  to  those  of  teleseismic  P.  A  brief  review  of  true  travel  time 
residuals  has  showed  however  that  the  scatter  of  residuals  is  much 
larger  for  teleseismic  P  than  for  Pn  and  Pg.  Since  true  residuals  do 
not  include  the  effect  of  source  mislocation,  these  values  reflect  true 
structural  heterogeneities  at  various  depths.  A  larger  scatter  in  P 
suggests  that  the  Earth  is  also  very  heterogeneous  in  the  lower  mantle. 
We  will  continue  to  investigate  this  in  the  future. 
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LOCATION  ERROR  ELLIPSES  FOR  RELS 


This  is  a  brief  note  which  attempts  to  clear  up  certain  recurring 
problems  in  our  assigning  confidence  limits  to  event  locations.  The 
problems  recur  because  we  haven't  adequately  documented  the  algorithm, 
so  there  is  always  uncertainty  and/or  faulty  memory  about  exactly  what 
it  does.  This  note  will  sketch  the  essentials  of  the  process  as  it 
stands  now,  so  that  the  necessary  changes  can  be  made  to  the  routine  in 
order  to  adjust  its  computation  to  those  which  are  required  by  RELS. 


For  an  initial  location  estimate  xQ>  the  arrival  time  (and  back 
azimuth)  residuals  are  expanded  as 

<$t  =  §  (x  -  xQ)  +  e 


(  A-l  ) 


where  we  assume  that 


cov(e)  =<ee’>  =  a2  T.  =  a2 


0  . 


(  A-  2  ) 


N 

ordinarily,  although  our  new  program  doesn't  require  that  E  be  diagonal. 
If  we  have  chosen  the  weight  w.  *  1/  a2,  to  be  the  true  variance  of  the 

t  h  li 

i  observation,  then  the  scale  factor  a2  =1.  The  leas t-sqaures 

solution  to  the  linearized  equation  is 

-*■  -»■  T  -1  -1  T  -1 

x  =  x  +  (B  Z  B)  B  Z  6k. 

0 . (A-*) 


We  see  that  back  azimuth  and  arrival  time  residual  may  both  appear  as 
components  of  6  even  though  they  don't  have  the  same  units,  since  they 
are  weighted  by  their  respective  £  priori  variances  so  that  all  the 
products  in  the  matrix  multiplication  have  units  of  km.  If  we  assume 
that  the  least-squares  estimate  x  has  a  mutivariate  Gaussian  distri¬ 
bution  about  the  true  location  x,  then 


(x  -  x) ' [o2(BT  I  1  5)  *1  1  (x  -  x)  <  x2q . 0 


( A-4 ) 


is  the  confidence  region  at  the  1-a  level  (we  choose  a  *  52),  or 


(x  -  x) '  (BT  T,  *  B)  (x  -  x)  _<  a2  \2p .  Q 


( A-  5 ) 


A- 2 


Here  q  is  the  number  of  dimensions  of  x.  For  the  conventional  location 
algorithm,  q  is  3  or  4,  depending  on  whether  it  is  run  in  the  depth-free 
or  depth-restrained  mode.  For  SIMUL,  however,  q  <  10,  since  we  now 
determine  coefficients  a„  ,  b_  ,  a_  ,  etc.,  which  describe  the  slopes 
and  intercepts  for  regional  phase  travel-time  relations  determined  Sl_ 
posteriori  (Chang  et  al.,  1981b).  If  we  do  not  believe  our  _a  priori 
estimate  a2  *  1,  we  replace  it  a  posteriori  with  the  calculated  variance 


“2  •  iM  t?t  -  ?<*  -  V1'  5*1  -  5  -  *o)1> 


(A-6) 


where,  for  the  final  iteration,  xQ  *  x  ,  so 


-1 

<5  £  1  <5 

t  -  t. 


( A-7) 


Note  that  we  still  assume  that  our  initial  assignation  of  the  relative 

weights  CT./CT.  is  correct;  obviously,  the  computed  location  x  as  well  as 
A  J 

the  confidence  region  will  change  if  we  decide  to  weight  P  the  same  as 

>s  n 

P,  for  example.  Now  a2  (N-q)/a2  has  the  y2  distribution  with  N-q 
degrees  of  freedom,  so  the  ratio  of  the  previous  X2statistic,  divided  by 

q  J 

q  degrees  of  freedom,  to  this  one,  divided  by  N-q,  has  the  F-distribu- 
tion  with  q,  N-q  degrees  of  freedom: 


(x  -  x)’[q2  (BT  r1  B)"1]"1  (X  -  x) '  1/c 
[a2  (N-q) /a2]/ (N-q) 


<  P 

—  q,N-q:a 


(A-3) 


and  the  confidence  region  becomes 


(x  -  x)'(BT  I-1  B)  (x  -  x)’  <  a  o2  F 


q ,N-q : a 


(A-9) 


The  confusion  arises  when  we  want  to  consider  spaces  of 
dimensionality  less  than  q.  For  an  m-dimensional  subspace  we  replace  x 
and  x  by  z  and  z,  their  projections  onto  the  subspace.  (Since  the 

components  of  x  are  x  “  [xyf  xx>  xfc ,  x^,  *apn . xbLg^  where  +  y 

north,  +  x  *  east,  ♦  t  *  future,  and  +  h  *  down,  we  will  conventionally 

perform  this  projection  into  3-  and  2-dimensional  space  by  simply 

dropping  the  final  7  then  8  components.  We  can,  of  course,  consider  4 

(or  higher)  -dimensional  confidence  regions,  as  was  done  m  the  initial 

version  of  LOCATION  by  Julian  (1973),  but  their  interpretation  is 

difficult  to  grasp,  and  it  will  be  better  to  consider  the  uncertainties 

*  (S  etc  .  in  terms  of  the  standard  errors,  as  we  shall  show  later.) 
t’  aPn 

Let  us  define 


-1  T  -1 
S  S  B  E  B 


(A-in) 


T  - 1  - 1 

so  that  S  *  (B  I  §)  is  a  qxq  matrix.  We  form  an  mxm  submatrix  E 
by  setting 

E..  *  S..  for  i  *  1 , . . . ,m  and  j  =  l,...,m. 

ij  lj 

Our  q-dimensional  chi-square  confidence  region  reduces  to 


(z  -  z)'  [a2  E]“  (z  -  z)  <  x^;01> 


( A-  1  1  ) 


Note  that  although  E  is  a  submatrix  of  S,  E  is  not  a  submatrix  of 
- 1  T  - 1 

S  =  B  £  B.  Watch  out  for  this  when  refering  to  Flinn's  classic 
article  on  error  ellipses  (1965),  which  has  problems  with  this 
distinction,  as  does  the  "corrections"  article  published  years  later. 

Now  no  changes  are  to  be  made  in  the  calculation  of  <52;  it  still  has  N-q 
degrees  of  freedom.  The  F-atataistic  confidence  region  in  m-dimensions 
is  thus 


fz  -  iy  E’1  (z-t)  »2 


( A-l 2) 


If  we  know  the  true  location,  as  in  the  case  for  NTS  shots,  we  can  now 
see  whether  this  confidence  region  around  the  computed  hypocenter  (m*3) 


A-4 


or  epicenter  (m*2)  does  in  fact  include  the  true  location.  Since  d2 is 

computed  by  the  program  and  is  stored  as  variable  "EMS"  in  common  block 

LOKCOM,  and  Since  S  is  passed  to  the  error  ellipse  subroutine  DMHLP 

via  the  argument  list,  we  can  simply  plug  the  appropriate  values  into 

this  formula  and  see  whether  the  left-hand  side  is  indeed  smaller  than 

the  right-hand  side.  (The  value  of  the  right-hand  sidee  is  computed 

automatically  by  DMHLP;  it  is  nno  longer  necessary  for  the  programmer  to 

input  the  appropriate  value  of  •)  We  can  also  now  determine  the 

relationship  between  the  2-d  error  ellipse  and  the  2-d  projection  of  the 

3-d  error  ellipsoid.  Although  E^^  is  a  submatrix  of  E  *2-d  *s 

not  a  submatrix  of  E  so  c^e  2-d  ellipse  differs  not  only  in  size 

(right-hand  side  of  the  inequality)  but  also  in  shape  (left-hand  side) 

from  the  2-d  projection  off  the  error  ellipsoid  (i.e.,  a  horizontal 

cross-section  of  the  vertical  cyclinder  containing  it).  [Note  that 

these  are  both  different  from  yet  a  third  2-d  "error  ellipse"  of  sorts, 

namely  the  intersection  of  the  3-d  ellipsoid  with  the  surface  plane  h  * 

0.)  It  is  thus  possible  that  the  true  hypocenter  lies  within  the  error 

ellipsoid  even  though  the  true  epicenter  lies  outside  the  error  ellipse, 

or  vice  versa.  As  an  example  of  the  first  situation,  we  note  that  if 

the  calculated  hypocenter  lies  directly  below  the  true  hypocenter  but  at 

a  much  greater  depth,  then  the  true  hypocenter  would  probably  lie 

outside  (namely,  above)  the  3-d  error  ellipsoid,  but  the  true  epicenter 

would  be  in  the  exact  center  of  the  2-d  error  ellipse.  As  an  example  of 

the  opposite  situation,  we  consider  the  special  case  in  which  one  of  the 

principal  (i.e.,  symmetry)  axes  of  the  3-d  error  ellipsoid  is  vertical. 

The  3x3  matrix  E,  .  may  then  be  written  as 
^  a 


( A- 1 3  ) 


on  e. 


so  we  have  E 


0  0 


( A-14) 


..H  . 


All  that  this  says  is  that  the  shape  of  the  2-d  projection  of  a 

vertically  symmetric  3-d  error  ellipsoid  is  the  same  as  that  of  the  2-d 

ellipse,  a  result  which  one  could  have  guessed  geometrically  and  which 

clearly  is  approximately  valid  even  if  the  3-d  ellipsoid  is  tilted 

slightly  from  the  vertical,  i.e.  ,  if  |  e13  |  <<  |  e33  |  and  |  e23  |  «  e33 

If  the  computed  hypocenter  and  the  true  hypocenter  are  at  the  same 

depth,  we  see  that  the  quantity  (Z-Z)'  E  1  (Z-Z)  is  the  same  for  the  3-d 

and  the  2-d  cases.  The  size  of  the  2-d  projection  of  the  3-d  ellipsoid 

differs  from  that  of  the  2-d  ellipse,  however;  in  fact,  the  ratio  of 

their  areas  is  given  by  3  F,  M  /2F.  =  1.3  for  N  .  In  this 

j , N-q  l  t N  —  q 

special  case,  then,  the  true  hypocenter  will  lie  within  the  3-d  error 

ellipsoid  but  the  true  epicenter  will  lie  outside  the  2-d  error  ellipse 

if  1.0  <(Z-Z)'  E_1  (Z-Z )  /  ( 2F-  62)<1.3. 

2,N-q 

We  can  now  compute  the  lengths  of  the  axes  of  the  m-d imens iona 1 
error  ellipsoid  and  the  spatial  orientations  of  these  axes.  In  order  to 
do  so,  we  rotate  from  the  fixed  y- ,  x- ,  and  h-axes  of  space  to  the 

symmetry  axes  of  the  error  ellipsoid.  In  this  new  corrdinate  system  the 
_ •  .-1  . _ 


matrix  E 


becomes 


(A- 13) 


where  X  ,  , 


X^  are  the  eigenvalues  of  E 


The  basis  vectors  of 


this  rotated  coordinate  system  are  of  course  the  eigenvectors  of  E  ; 
since  they  are  orthonormal,  their  components  (expressed  in  the  y-,  x- , 
h-space)  are  their  direction  cosines  with  respect  to  the  fixed  axes.  Tn 
the  rotated  coordinate  system,  which  we  denote  as  Sj ,  C2>  ,  the  3-d 

error  ellipsoid  is 


+^2  +  X3  ^3  -  F3’N-q:  a. 


(A- It 


with  a  similar  equation  for  the  2-d  error  ellipse,  where  we  have  now 
simplified  matters  by  setting  the  origin  of  the  coordinate  system  at  the 


A- 6 


computed  location  Z.  The  length  of  the  i^-axis  of  the  ellipsoid  is 
given  by  setting  and  equal  to  zero  in  this  expression: 


sw 


3,N-q 


A. 

a  1 


(  A-  l  7  ) 


and  similarly  for  the  lengths  of  the  t;  -  and  C^-axes.  We  can  thus  find 
the  lengths  and  orientations  of  the  error  ellipsoid  axes  by 
diagonalizing  E  our  program  uses  the  subroutine  EIGEN  to  do  this. 

The  SOL  report  on  HYPO  (Flinn,  1965b)  indicates  that  a  geometrical 
manipulation  is  used  in  that  program  to  accomplish  this,  but 
diagonalizing  E  *  seems  to  be  more  straightforward,  especially  in  the 
3-d  case. 


Even  though  we  now  know  the  lengths  of  the  ellipsoid  axes,  we  still 
do  not  know  the  maximum  range  (at  the  1-  a  confidence  level)  in  each 
direction  of  the  possible  values  of  x  .  For  example,  in  the  diagram 
below,  we  know  the  lengths  of  the  semi-major  and  semi-minor  axes  and 
^  and  their  azimuths  S  ^  and  but  we  have  not  yet  determined  ^  y  and 

6  x •  ( Figure  1 ) ; 


V 


For  the  moment  let  us  scale  the  ellipse  to  a  unit  size  by  replacing  the 
right-hand  side  of  its  equation  by  unity: 


(z  -  z)  E  1  (z  -  z)  *  1 


(  a- i a ) 


which  we  expand  as 


'11 


y2  +  2b 1 2  *y  +  b 2 2  x" 


1 


A- 7 


( A—  1  9  ) 


where  E 


Differentiating  with  respect  to  x, 
we  obtain 


2  bn  yy'  +  2  bJ2  y  +  2  b12  xy'  +  2  h  22  x  =  0 


(A-20) 


where  y'  =  —  .  At  the  maximum  extension  in  the  y-direction  of  the 

dx 

ellipse  5y,  we  have  y'r  =  0,  so 


2  b12  S'  «y  *  2  b22  X  fiy  *  0 


‘  (-b12,b22>>’fv. 


(A-211 


Substituting  into  the  equation  for  the  ellipse,  we  have 


bll  V6v  _2(b12/b22)  y6v  +  (b12/b22)  V6v  =  1 


( A-2  2 ) 


(bll  ”  b12  /b22^  Vfv  “  1 


v(,v  =  ~  bl2  ^b22  . 


( A-2  3 ) 


The  maximum  northward  extent  of  the  error  ellipse  is  thus  given  in  terms 

-1  ,  r 

of  E  ;  a  similar  analysis  can  be  used  to  determine  °x  or,  in  the  3-d 
case ,  & h . 


T  - 1  - 1 

A  much  faster  approach  involves  using  the  matrix  5  =  ( B  j  B) 
rather  than  the  matrix  g  '  To  see  this,  consider  the  case  in  which  the 
submatrix  E  is  in  fact  the  entire  matrix  S.  (This  would  happen  if  the 
location  were  performed  with  all  variables  restrained  except  for  x  and 


y  .  )  I  f  S 


,  then 


E~l  =  S_1 


all  a22  ~  a12  l  -a 


22 

f 

CN 
r— 4 

CO 

1 

( bn 

b12 

12 

a22  ] 

\b12 

b22 

(A-24) 


A-8 


2 

6v  *  1/ *4)^2  -  b12/b  22 

=  A12  //h1]  b22  -  bf2 
_1 _ 

=  /au  /[*^n  a22  "  a12  ^11  a22  "  a12^ 


( A-25) 


In  general 


=  {[cov  (x))^}  ^ 

1  n  \  1  ml.  _  _ -  - 


( A- 26 ) 


where  cov(x)  =  a2  (B  £  B)  .  The  maximum  extents  in  each  direction 
of  the  unit  error  ellipse  (the  "standard  errors")  are  thus  simply  the 
square  roots  of  the  diagonal  elements  of  S.  We  must  now  scale  these 
values  by  replacing  unity  on  the  right-hand  side  of  the  defining 
equations  by  some  constant  K;  the  standard  error  will  then  be  multiplied 
by  the  factor  /"kT  In  our  analyses  of  the  2-d  and  3-d  confidence 
regions  we  have  used  K  =  2  3%,,  „  and  K  =  3d2  F,  „  .  As  we  have 

pointed  out,  this  depends  on  the  dimensionality  of  the  confidence 
region;  for  the  3-d  error  ellipse,  the  maximum  ranges  6  x  and  <5  y  were 
shown  to  be  times  larger  than  <5x  and  <$  y  for  the  2-d  error  ellipse 

in  a  special  case.  The  best  way  of  treating  the  uncertainty  in  a  single 
variable  is  to  regard  it  simply  as  a  one-dimensional  confidence 
ellipsoid.  If  we  assume  that  the  variance  is  known  a  priori ,  we  thus 


<5v  -  {[cov  (x)]..  }  •  /~^2  y2 

11  1:  a 


T  -1  -1  1/2 

{[(B1  7  B)  }  •  1.960. 


(4-27) 


(if  a2  *  1,  which  we  define  it  to.)  This  is  of  course  just  what  you 
would  expect  for  the  95%  confidence  region,  namely  that  it  is  1 . 96 
standard  deviations.  Note  that  this  a  priori  estimate  depends  solely  on 
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the  geometry  of  the  stations  and  the  hypocenter,  and  it  does  not  depend 
on  the  measured  arrival  times.  We  are  therefore  more  interested  in  the 
a  pos  ter iori  estimate 


A  1/2 


6y  =  {[cov  (x)]  }  •  fo* 

11  1 ,  N-q  :  a 


T  -1  -1  1^2  _ 

=  C[B  I  B)  ]  }  8  fF~~ 

~  11  l,N-q:a. 


(  A—  2  8  ) 


This  is  the  best  way  to  handle  the  uncer  taint  ies  5  h  and  <5 1 .  The  program 
currently  calculates  the  standard  errors  and  the  variances  d2,  so  it  can 
trivially  be  adjusted  to  compute  &  h  and  6t  in  this  manner.  We  propose 
that  the  output  of  the  RELS  lcoation  error  subroutine  be: 

the  a_  posteriori  variance 

the  lengths  of  the  semi-major  and 
semi-minor  axes  of  the  confidence 
e 1 1  ipse 

the  azimuths  (in  radians  clockwise 
from  north?)  of  these  axes 


62 

?l(2-d),  ?2(2-d) 
?i(2-d),  C2(2-d) 


£h,  6t  the  standard  errors  in  depth  and 

origin  time,  multiplied  by  the 
appropriate  factors. 


The  parameters  q(2_d),  C2(2-d)’  5h(l-d)  and  lSt(l-d)  C3n  be  evaluated  for  *2 

and/or  F-statistic  cases.  Optional  parameters  to  return  are  C^d-d)' 

^  (3-d )  *  ^3  (3-d)  *  A  =  *  *  %(2-d)  ^2  ( 2-d )  ’  and  V  =  ^  l(3-d)  *  c2(3-d)'  ?3(3-d) 

in  the  y2 and/or  F-statistic  cases.  We  recommend  that  the  angles  £  be 

expressed  relative  to  the  fixed  y-,  x- ,  and  h-  axes  rather  than  as  Euler 

angles  within  the  ellipsoids,  as  was  done  by  Julian  (1973). 
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